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Thermodynamics of the unitary Fermi gas at finite temperature is investigated from the per- 
spective of the expansion over e = 4 — d with d being the dimensionality of space. We show that 
the thermodynamics is dominated by bosonic excitations in the low temperature region T -C T c - 
Analytic formulas for the thermodynamic functions as functions of the temperature are derived to 
the lowest order in e in this region. In the high temperature region where T ~ T c , bosonic and 
fermionic quasiparticles are excited. We determine the critical temperature T c of the superfluid 
phase transition and the thermodynamic functions around T c to the leading and next-to-leading 
orders in e. 
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I. INTRODUCTION 

The Fermi gas with zero-range interaction at infinite 
scattering length [1-3] , frequently referred to as the uni- 
tary Fermi gas, has attracted intense attention across 
many subfields of physics. Experimentally, the system 
can be realized in atomic traps using the Feshbach reso- 
nance and has been extensively studied [4-10]. Since the 
fcrmion density is the only dimensionful scale of the uni- 
tary Fermi gas, its properties are universal, i.e., indepen- 
dent of details of the interparticle interaction. The uni- 
tary Fermi gas is an idealization of dilute nuclear matter 
and may be relevant to the physics of neutron stars [11]. 
It has been also suggested that its understanding may be 
important for the high-T c superconductivity [12]. 

The austere simplicity of the unitary Fermi gas implies 
great difficulties for theoretical treatments, because there 
seems to be no parameter for a perturbation theory. Re- 
cently, we have proposed a new approach for the unitary 
Fermi gas based on the systematic expansion in terms 
of the dimensionality of space [13, 14], utilizing the spe- 
cialty of four or two spatial dimensions in the unitarity 
limit [15]. In this approach, one would extend the prob- 
lem to arbitrary spatial dimensions d and consider d is 
close to four or two. Then one performs all calculations 
treating £=4-dore = d- 2asa small parameter for 
the perturbative expansion. Results for the physical case 
of three spatial dimensions are obtained by extrapolating 
the series expansions to e (e) = 1, or more appropriately, 
by matching the two series expansions. 

We used this e expansion around four spatial dimen- 
sions to calculate thermodynamic functions and fermion 
quasiparticle spectrum in the unitarity limit and found 
results which are quite consistent with those obtained by 
Monte Carlo simulations and experiments [13, 14]. This 
e expansion has been successfully applied to atom-dimcr 
and dimcr-dimer scatterings in vacuum [16]. Thus there 



are compelling reasons to hope that the limit d — > 4 is not 
only theoretically interesting but also practically useful, 
despite the fact that the expansion parameter e is one at 
d = 3. Very recently, the phase structure of polarized 
Fermi gas near the unitarity point has been investigated 
based on the e expansion [14, 17]. 

In this paper, we extend our approach to investigate 
the thermodynamics of unitary Fermi gas at finite tem- 
perature T, below and above the critical temperature 
T = T c of the superfluid phase transition. So far, such 
a problem has been studied relying on the mean-field 
description with fluctuations [18-28], the virial expan- 
sion [29-31], or the Monte Carlo simulations [32-36]. The 
main advantage of our approach is that e = 4 — d provides 
the small parameter of the perturbative expansion, and 
thus, analytic and systematic study below and above T c 
is possible. We also note that the critical dimension of 
a superfluid-normal phase transition is also four, which 
makes weak-coupling calculations reliable at any temper- 
ature for the small e. 

First, we review the e expansion for the unitary Fermi 
gas around four spatial dimensions in Sec. II. Then we 
discuss that the thermodynamics in the low tempera- 
ture region T <C T c is dominated by bosonic excitations. 
Analytic formulas for the thermodynamic functions as 
functions of the temperature to the lowest order in e are 
shown in Sec. III. The behavior of the thermodynamic 
functions above T c to the leading and next-to-leading or- 
ders in e is discussed in Sec. IV. In particular, we put em- 
phasis on the determination of the critical temperature T c 
by matching the e expansion with the expansion around 
two spatial dimensions. The thermodynamic functions at 
T c are also calculated, which are compared to the results 
from recent Monte Carlo simulations. Finally, summary 
and concluding remarks are given in Sec. V. 
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II. REVIEW OF e EXPANSION 

Here we briefly review the e expansion for the unitary 
Fermi gas around four spatial dimensions and quote some 
results at zero temperature to the leading order in e just 
for the convenience of later discussions at finite temper- 
ature. The detailed account of the e expansion is found 
in Refs. [13, 14]. 



A. Lagrangian and Feynman rules 

The extension to finite temperature T follows from the 
prescription of the imaginary time formalism. The sys- 
tem under consideration is described by the sum of fol- 
lowing Lagrangian densities (here and below h = 1 and 
ks = 1): 
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The propagators of fermion and boson are generated by 
Cq. The fermion propagator G is a 2 x 2 matrix, 
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where e p = p 2 /2m is the kinetic energy of nonrelativis- 
tic particles and </>o chosen to be real is the condensate in 

the superfluid ground state. E p = 



■q is the ex- 



citation energy of the fermion quasiparticle. The boson 
propagator D is given by 
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co n — 2irT(n + ^) and v n — 2irTn are discrete Matsubara 
frequencies for fermion and boson with an integer n — 
0,±1,±2,---. 

The unitary Fermi gas around four spatial dimensions 
is described by the weakly interacting system of fermionic 
and bosonic quasiparticles, whose coupling g ~ e 1 ' 2 in C\ 
is given by 
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FIG. 1: Restoration of naive e counting for the boson self- 
energy. The fermion loop in (c) goes around clockwise and 
counterclockwise. Solid (dotted) lines represent the fermion 
(boson) propagator — G (-D), while the cross in (c) represents 
the (j, insertion to the fermion propagator. 



The chemical potential fi ~ e(f>o in C\ is treated as a 
small perturbation in our formulation. We define the 
boson chemical potential as 



Mb = 2/x 
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In the dimensional regularization we use, when Co is neg- 
ative, —g 2 /co ~ £b gives the binding energy of the boson 
to the leading order in e. Throughout this paper, we con- 
sider the vicinity of the unitary point where £b ~ e<^o- 

Finally the two additional vertices for the boson prop- 
agator in Li play a role of counter terms so as to avoid 
double counting of certain types of diagrams which are 
already taken into Cq and C\. The first vertex in the 
momentum space is given by 



MPo,p)=Po~f- (8) 



B. Power counting rule of e 

The power counting rule of e is summarized as follows. 

1. We consider ^/0o ~ e and regard 4>o as 0(1). 

2. For any Green's function, we write down all Feyn- 
man diagrams according to the Feynman rules us- 
ing the propagators from Cq and the vertices from 
Ci. 

3. If there is any subdiagram of the type in Fig. 1(a) 
or Fig. 1(c), we add the same Feynman diagram 
where the subdiagram is replaced by a vertex from 
C 2 , Fig. 1(b) or Fig. 1(d). 



4. The power of e for the given Feynman diagram will 
be 0(e Na / 2+N >*), where N g is the number of cou- 
plings g and is the number of chemical potential 
insertions. 
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5. The only exception is the one-loop vacuum diagram 
of fermion with one fj, insertion, which is O(l) in- 
stead of the naive O(e). 

This power counting rule holds in the low temperature 
region where the condensate is still large compared to 
the chemical potential m/^o ~ e, while it breaks down 
near the critical temperature because <fo — > at T — > 
T c . In the high temperature region T ~ T c , a minor 
modification of the power counting rule is necessary as 
we discuss at the beginning of Sec. IV. 
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C. Leading order results at T — 

Here we show only lowest order results on the thermo- 
dynamic functions at zero temperature, which is equiva- 
lent to the mean-field approximation. The results up to 
the next-to-leading order in e are found in Refs. [13, 14]. 
The effective potential to the leading order in e is 
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from which the condensate follows as 
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Note that the previously made assumption fJ,/<po = 0(e) 
is now checked. Then the pressure at zero temperature 
is given by 
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From the fermion number density 
N 
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we define the Fermi energy through the relationship in 
the ideal Fermi gas in d spatial dimensions; 



FIG. 2: Boson's self-energy diagrams contributing to the or- 
der O(e). The vertex IIo from C2 needs to be added to the 
second diagram. The last diagram gives the off-diagonal part 
of the self-energy. Other elements are given by Il22(p) = 
IIii(— p) and Il2i(p) = Uuip)*. 



III. THERMODYNAMICS BELOW T c 

Now we investigate the thermodynamics of the Fermi 
gas at finite temperature near the unitarity limit. At zero 
temperature, we found that there exist two difference en- 
ergy scales in the system; the scale of the condensate 4>o 
and that of the chemical potential \x ~ ecbp <C <j>o. Ac- 
cordingly, we can consider two temperature regions where 
the unitary Fermi gas exhibits different thermodynamics. 

One is the low temperature region where T ~ eahp. In 
this region, the energy gap of the fermion quasiparticle 
excitation A ~ cbp is still large compared to the tem- 
perature. Therefore, thermal excitations of the fermion 
quasiparticle are exponentially suppressed by a factor 
e -A/T ^ e ~ 1 / e . The thermodynamics in this region is 
dominated by the bosonic phonon excitations. The other 
temperature region is the high temperature region where 
T ~ </>o- (0o represents the condensate at zero tempera- 
ture.) In this region, the condensate decreases and even- 
tually vanishes at the critical temperature T c . Fermions 
and bosons are equally excited here. We defer our dis- 
cussion on the high temperature region to Sec. IV and 
concentrate on the thermodynamics in the low tempera- 
ture region T -C T c in this section. 
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Then the energy density is given by 



E Q = f iN - P = 
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where Sb = —g 2 /co > is the binding energy of the 
boson. The chemical potential at the fixed Fermi energy 
£f is obtained from Eqs. (10) and (13) as 
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A. Phonon spectrum 

The thermodynamics in the low temperature region 
T -c cf>Q is dominated by the phonon excitations. In order 
to determine the phonon spectrum, we first study the bo- 
son self-energy at zero temperature. To the order of 0(e) , 
there are three types of contributions to the boson self- 
energy as depicted in Fig. 2. In addition to the chemical 
potential insertion /jb — 2m + £b, the one-loop diagrams 
contribute to the diagonal part iln(p) and off-diagonal 
part IIi2(p) of the boson self-energy. The vertex Ho(p) 
from £2 is necessary for Hn(p) according to the power 
counting rule described in Sec. II. Other elements can be 



4 



obtained from II 22 (p) = IIn (— p) and n 2 i(p) = IIi 2 (p)*. where U (p) is denned in Eq. (8) and n a (p) is given by 

The diagonal part of the boson self-energy of 0(e) is 
given by the sum 

n 11 (p) = n ( P ) + n a (p), (16) 



(E k _E + e k _R){E k+ i 
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Here we introduced the shorthand notation, 

dk 
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(18) 



Since we are interested in physics at the scale of temper- 
ature T <C <fro, it is sufficient to evaluate the self-energy 
when the external momentum is small compared to the 
condensate (j>o,£ P ) ~ T -C </>o- Expanding IIii(p) in 
terms of p/</>o and performing the fc integration with the 
use of the formula 



T(a)T{fi-a) 



(19) 



we obtain 
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= -e0o + O(e 2 ). (20) 



Similarly, the off-diagonal part of the boson-self energy 
of 0(e) is given by 
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Expanding II 12 (p) in terms of p/<f>o and performing the 
integration over fc, we obtain 



ni 2 (0) =g : 



0(e 2 ). 



(22) 



As a result of the resummation of these self-energies, 
the resummed boson propagator V is expressed by the 
following 2x2 matrix: 



2?(Po,p) = 

-n 2 i 



-n 12 

^(-p)- 1 + MB - n 22 



(23) 



where II 22 (p) = IIii(-p) and n 2 i(p) = ITi 2 (p)*. Using 
the self-energies calculated above, the dispersion relation 
of the boson w p h(p) can be obtained by solving the equa- 
tion det[.D _1 (w,p)] = in terms of u as 



^Ph(p) = y (y - PB + e^o) ( y - MB + 2e0 o ) ■ (24) 

Note that this expression is valid as long as e p -C 4>o be- 
cause of the expansions made to evaluate the boson self- 
energies. Substituting the solution of the gap equation 
at zero temperature in Eq. (10), Mb = e^o, the phonon 
spectrum is determined to be 



= Vt (t 



£00 



(25) 



For the small momentum e p <C e</>o, the dispersion 
relation becomes linear in the momentum as cj p h ~ 
c s |p|, remaining gapless in accordance with the Nambu- 
Goldstone theorem. The sound velocity of phonon c s is 
given by 



£T (26) 




where vf — (2£p/m) 1//2 is the Fermi velocity. For the 
large momentum e p 3> e</> , the dispersion relation ap- 
proaches that of the free boson as w p h — £ P /2. 



B. Effective potential and condensate 

At finite temperature, the phonon excitations con- 
tribute to the effective potential, and consequently, the 
magnitude of the condensate decreases. The tempera- 
ture dependent part of the effective potential Vr((f>a) to 
the lowest order in e is given by the one-loop diagram 
of the boson with the resummed propagator in Eq. (23) 
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FIG. 3: One- loop diagram of boson contributing to the effec- 
tive potential at finite temperature. The dotted double line 
represents the resummed boson propagator T> in Eq. (23). 



where V e ff(<po) is the effective potential at zero tempera- 
ture in Eq. (9). To the leading order in e, tfo at T <C e^o 
is given by 



while at e^o CTc^o, 



8C(3)T 3 



n 2 T 2 



(33) 



(34) 



(Fig. 3): 
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The condensate in total is 4> — + 4>t, which de- 
creases as the temperature increases. Note that since 
(j)T <C e(j) , the leading part of the condensate docs 
not change in the temperature region considered here 
T <C 0o- The effective potential is given by the sum 
of the zero temperature and finite temperature parts; 

V e ffW>0 + <h) + v t(0o + <h) - VcsiM + v T (M- 



The n-point interaction vertex among phonons (p n is of 
the order e ra//2 and appears in the effective potential only 
at higher orders. Then the contribution of Vt(<Po) to the 
gap equation is 



dV T (M 
d(j) 



/B(Wph) 



(28) 



where /b(x) = l/(e x / T — 1) is the Bose distribution func- 
tion and du> p h/d(j)o is given from Eq. (24) by 



duj ph _ 3e p + 2e0 o 



4w 



P h 



(29) 



There are two limiting cases where the integration over 
p in Eq. (28) can be analytically performed. Since the 
integral is dominated by the integration region where 
e p ~ T, we can approximate the phonon spectrum by 
its linear branch u p h(p) — c s |p| when the temperature is 
very low T <c e(f>Q. In this case, the integration over p at 
d = 4 leads to 



dVr 



6o) ^ 8C(3)T /mT 

I 4>0 \ 27T 
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On the other hand, when the temperature is located in 
the intermediate region e(f>Q C T < ^, the phonon 
spectrum can be approximated by its quadratic branch 
Wph — Sp/2. In this case, the integration over p in 
Eq. (28) at d = 4 results in 



dv T {<M ^ c KH 2 



(31) 



Now, from the gap equation d[V e g((p) + Vr(4>)]/d<j> = 
with <j) = 4>q + (f>T, one finds the temperature dependent 
correction of the condensate 4>t satisfies 



dV cS (fo)\ , dV T ( 

VT H 7T, 
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C. Thermodynamic functions at low temperature 

The temperature dependent part of the pressure P p h 
in the low temperature region T <C <fro is given from the 
effective potential in Eq. (27) by 



Pph = -V T (M = -T j 
Jv 



In 



_ e -w pi (p)/T 



(35) 



The phonon contributions to the fermion number den- 
sity, the entropy density, and the energy density are 
respectively computed from the thermodynamic rela- 
tions, N ph = dP ph /dn, S ph = dPp h /dT, and E ph = 
M^Vph + TSph — P p h- Here we show analytic expressions 
for these thermodynamic functions in the two cases where 
the analytic evaluation of the p integration in Eq. (35) is 
available. 

When the temperature is very low T <C e</>o, only the 
linear branch of the phonon spectrum oj p h(p) — c s \p\ 
is important to the thermodynamic functions. In this 
case, the integration over p at d — 4 can be performed 
analytically to lead to 



12tt 2 C(5) T 5 _ 12C(5) m 2 T 5 

(27T)4 C 4 " n 2 ( £ O )2' 
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Note that the exponent of T in general spatial dimension 
d is Pph ~ T d+1 . Accordingly, we obtain the phonon 
contributions to the fermion number density, the entropy 
density, and the energy density: 
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Since actual experiments or simulations are performed 
with the fixed fermion density, it is useful to show the 
thermodynamic functions at fixed N instead of fixed fx. 
From Eqs. (10), (12), and (37), we find the chemical po- 
tential for the fixed fermion density increases as a func- 
tion of the temperature as 



2^5 

Ai = /i + 48C(5)— 4, 



(40) 



where Mo represents the chemical potential at zero tem- 
perature in Eq. (15). Normalizing zz by the Fermi energy 
in Eq. (13), we have 



H _ Mo , 3C(5) /2TV 
£ F e F 2e 3 \e F J 
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The other thermodynamic functions for the fixed fermion 
number density are given by 
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These expressions are valid in the low temperature region 
where T -C efo and hence T/e F -C e 3 / 2 . The exponent 
of T is different from that in three spatial dimensions be- 
cause we are expanding around four spatial dimensions. 
The correct exponent at d = 3 is recovered if we resum 
logarithmic corrections ~ (— elnT)™ to infinite orders, 
which arc formally higher orders in e and are not shown 
in the formulas above. 

On the other hand, when the temperature is located in 
the intermediate region e<j)Q <C T <C 4>q, we can expand 
the phonon spectrum w p h(f>) in terms of e<j>o/s p up to its 
first order 



^ph(p) - y + — • 



(45) 



In this case, the integration over p in Eq. (35) at d = 4 
can be performed analytically again to result in 
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Note that the exponent of T in general spatial dimen- 
sion d is P ph ~ T d/2+1 + T d/2 . Accordingly, we obtain 
the temperature dependent parts of the fermion number 
density, the entropy density, and the energy density: 
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From Eqs. (10), (12), and (47), we find the chemical 
potential for the fixed fermion density increases as a func- 
tion of the temperature as 

7T 2 T 2 

Li = lM) + e 2 — — . (50) 
6 0o 

Normalizing fj, by the Fermi energy in Eq. (13), we have 

(51) 
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The other thermodynamic functions for the fixed fermion 
number density are given by 
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These expressions are valid in the intermediate temper- 
ature region where e</> <C T C 4 and hence e 3 / 2 -C 
T/e F -C e 1 ' 2 . The exponent of T is different from that 
in three spatial dimensions because we are expanding 
around four spatial dimensions. The correct exponent 
at d — 3 is recovered if we resum logarithmic corrections 
<~ (— | In T) to infinite orders, which are formally higher 
orders in e and are not shown in the formulas above. 



D. Effective potential near T c 

At the end of this section, we study the behavior of 
the condensate as a function of the temperature for 
a given fi. The effective potential to the leading order 
in e is given by one-loop diagrams of fermion with and 
without one /x insertion, which is equivalent to the mean- 
field approximation. Since the critical dimension of the 
superfluid-normal phase transition is four, the mean-field 
approximation remains as a leading part at any temper- 
ature in the limit d — > 4. Then the leading contribution 
to the effective potential at finite temperature is given by 

VeS (<t>) = -%<?- / Up ~ 1T» + 2Tln (l+e- B " /T ) 
9 Jp L 

(55) 

where f F (x) = l/(e x ' T + 1) is the Fermi distribution 
function. 

For the low temperature T -C (j), we can neglect the ex- 
ponentially small factor e~ E vl T <c 1 and the integration 
over p reproduces the effective potential at zero temper- 
ature in Eq. (9). In the opposite limit where 4> <§; T, we 
can expand V e s((j>) in terms of <p/T to lead to 

V eS (<f>) ~ Kff(0) (56) 
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FIG. 4: The behavior of the condensate <f> as a function of 
the temperature T to the lowest order in e. 0o = Mb/e is the 
value of the condensate at T — and the critical temperature 
is located at T c /cj> = 1/(2 In 2) « 0.721348. 



IV. THERMODYNAMICS ABOVE T c 
A. Power counting rule of e near T c 

In the e expansion at zero or low temperature, the 
chemical potential is small compared to the condensate 
(i ~ e<po and we made expansions in terms of fi/4>o ~ e 
as well as the small coupling g ~ e 1 / 2 . Near the critical 
temperature, the ratio fj,/(f> is no longer small because the 
condensate <j> vanishes at T = T c , while n/T c is O(e) as 
it is clear from T c = /i B /(e21n2) in Eq. (57). Therefore, 
we can still treat the chemical potential as a small per- 
turbation near T c and the same power counting rule of e 
described in Sec. II B holds even above T c just by replac- 
ing (f> with T. Hereafter we regard T ~ T c as 0(1). 



where = 2^+£b- From the coefficient of the quadratic 
term in <j>, we can read the critical temperature T c of the 
second order phase transition to the leading order in e as 

and the value of the condensate (j) just below T c as 

<j) 2 = 8T(T C -T)ln2 + 0(e). (58) 

The condensate <j> in the intermediate range of the 
temperature is obtained by solving the gap equation 
dV cS /d<f> = 0: 

<t>-^+ Jde p ^ME p ) = 0. (59) 

The numerical solution of the gap equation as a function 
of T is shown in Fig. 4. 

Finally we note that critical exponents for an 0(2)- 
symmetric theory will be in principle recovered if we re- 
sum logarithmic corrections ~ (e m In \T — T c \) n appear- 
ing at higher orders due to the infrared physics of bosons 
with a zero Matsubara frequency. [We can see such a log- 
arithmic correction, e.g., in Eq. (69).] Resumming these 
logarithmic corrections to infinite orders, the critical ex- 
ponent shifts from its mean-field value by \T — T c \# e to 
produce the known results for the XY universality class 
in the usual e expansion. In other words, we cannot see 
the shift of the critical exponent at finite orders in our e 
expansion. 



B. Boson's thermal mass 

First we study the self-energy of boson at T > T c . The 
leading contribution to the self-energy is the chemical 
potential insertion /xb as well as the one-loop diagram 
IIii shown in Fig. 2: 

Uu(iv,p) - Yl (iv,p) 

n Jk 

_ _ 2 f 1 ~ Mffc-g/g) ~ h{ek+ P /2) , . 

For the zero Matsubara frequency mode v n — at the 
small momentum e p -C T, we have 

nn(0,0) =g 2 f Mf^l = er2 ln2 + (9(e 2 ). (61) 
Jk £fe 

Therefore, the zero Matsubara frequency mode has the 
non-negative thermal mass Hr = eT21n2 — ^b ~ O(e) 
at T > T c . The condition of the vanishing thermal mass 
Ut = gives the critical temperature T c = ^ B /( e 2 ln2) + 
0(e) equivalent to Eq. (57). As we will see below, at a 
sufficiently high order in the perturbation theory near 
T c [0(e 2 ) or 0(e) compared to the leading term in the 
pressure or fermion density], a resummation of the bo- 
son self-energy is needed to avoid infrared singularities 
appearing in the zero Matsubara frequency mode. 



C. Pressure 



Now we calculate the thermodynamic functions at T > T c to the leading and next-to-leading orders in e. There 
are three types of diagrams contributing to the pressure up to the ncxt-to-leading order in e as depicted in Fig. 5; 
one-loop diagrams of fermion (boson) with and without one \x (i^b) insertion and two- loop diagram with a boson 
exchange. Note that at T > T c , the boson's one-loop diagram contributes as O(l) as well as the fermion's one-loop 
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+ t 



+ 



0(1) 




FIG. 5: Three types of diagrams contributing to the pressure up to the next-to-leading order in e. Each /i Qi/b) insertion to 
the fermion (boson) line reduces the power of e by one. 



diagram. Then the pressure from the one-loop diagrams is given by 



Pi 



2Tln (l + e - £ »/ T ) - Tin (l - e -^/ 2T ) + 2 M / F (e p ) + Mb/b(£ p /2) 



T 



11 91n2C(3) + llC'(3) 7T 2 M , 2tt 2 M b 



d/2 



The contribution from the two-loop diagram to the pressure, which is 0(e), is given by 

Pi = g 2 T 2 ^2 / Gii(iuj n ,p)G 2 2(i^m,q)D(iuj n - iu m ,p - q) 
^ f fF(e P )fF(e q ) + [f F (e p ) + Me q )}f B (s p - q /2) 



l P q £j> "F £ P —q/2 

The numerical integrations over p and q result in 



d/2 

P2 = -C P e [^-) T, 



(62) 



(63) 



(64) 



where Cp « 8.4144 is a numerical constant. From Eqs. (62) and (64), we obtain the pressure up to the next-to-leading 
order in e as 



P = P 1 +P 2 =T 



11 91n2C(3) + 11C'(3) 7T 2 /x 2tt 2 Mb 

yC(3) s e -^ + It + TY 



fmT\ 



d/2 



(65) 



The entropy density S and the energy density E to the same order can be computed from the thermodynamic relations 
S = dP/dT and E = fiN + TS - P. 



D. Fermion number density 



The fermion number density to the next-to-leading order in e cannot be obtained simply by differentiating the 
pressure in Eq. (65) with respect to the chemical potential N — dP/d/i. Since the pressure to the leading order in 
e does not depend on fi and the /i derivative d/dj-i ~ 1/e enhances the power of e by one, we need to compute the 
one-loop diagrams with two /i (pb) insertions and the two-loop diagrams with one /i (pb) insertion. The fermion 
density from the fermion's one-loop diagrams is given by 



N F 



V 

^ 7T 2 hl2 

~6 



6C'(2) 



2 In 2 



-li 



(mT\ 



d/2 



(66) 



12 T 

On the other hand, the boson's one-loop diagrams contribute to the fermion density as 



N B = 2 



h(e P /2) - ^/ B ( £p /2)/ B (- £p /2) 



(67) 
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n 



FIG. 6: A type of diagrams resummed to resolve the infrared singularity in Eq. (67). 



Apparently, the last term has an infrared singularity because the Bose distribution function behaves as f B (e p /2) ~ 
2T/e p at small momentum e p <C T. In order to resolve this infrared singularity, the resummation of the boson 
self-energy Ht = Tin — Mb ~ 0(e) is needed at the small momentum region £ p <~ fi. After resumming a type of 
diagrams shown in Fig. 6, we can rewrite N B within the accuracy we are working as 



2 /"[/ B (e P /2 + n T )-e21n2 



/ B (£ P /2)/ B (-£ p /2)] + 0(e 2 ), 



(68) 



where we used Tln/T — e21n2. Now the first term is infrared finite, where the boson's thermal mass Tlx plays a role 
of an infrared cutoff. Integrating over p and expanding up to the next-to-leading order in e, we have 

d/2 



2 f f B (e p /2 + n T ) 



4^ 
~3~ 



The logarithmic term ~ In Tlx /T appears as a consequence of the resummation. The second term in Eq. 
is still infrared divergent, will cancel with the infrared singularity existing in the two- loop diagram. 
The contribution from the two-loop diagrams to the fermion density is given by 



N 2 = -g* 



pk 



/F(£fc+p/2)/F(-£fc+p/2) [Mffc-g/g) + /b(£ p /2)] - 2/ F (£ fc+ p/2)/ B (£p/2)/B(-£p/2) 



(69) 
which 

(70) 



The second term in the numerator contains the infrared singularity at small e p . Extracting the divergent part, we 
can rewrite N 2 as 



N 2 = -g 2 



M£fc+p/2)/F(-£fc+p/2) [fc(e k _ p /2) + / B (£ P /2)] - 2 [/ F (£fc+ P /2) - Mefc)] /b(£ p /2)/b(-£ p /2) 



pk 

2.9 2 



/F(£fc) 



/b(£ p /2)/b(-£ p /2). 



(71) 



l p k £ kT 

Now the first term is infrared finite. One finds the k integration in the second term can be performed to lead 
to n n = eT21n2 in Eq. (61), which exactly cancels out the infrared divergent part in Eq. (68). The numerical 
integrations over k and p in the first term result in 

d/2 

N 2 = ~C N e +e41n2 / / B (£p/2)/ B (-£ P /2), (72) 



fmT\ 1 f 
-C N e ^_ j +e41n2 J / B ( £p /2)/ B (-£ p /2), 



where Cjv ~ 1.92181 is a numerical constant. 

Gathering up all contributions, Eqs. (66), (68), and (72), the fermion number density to the leading and next-to- 
leading orders in e is given by 



N = N F + N B + N 2 

"3tt 2 37r 2 ln2 + 18C'(2) 



2in2 / , n T \ n T 

eCjv + ^M-Sf 1-ln — J — 



2 4 

We define the Fermi energy £p through the relationship in Eq. (13) as 

3C'(2) 



£p 

T 



3tt 2 



1 



27 - 3 - 2 In 2 



2tt 2 



3^ 



3tt 2 



2 In 2 n 
3tt 2 T 



s 



3^ 2 K 1 -^— ~ 



(73) 



(74) 



The logarithmic correction (TTt/T) \iiTTt/T ~ elne is a consequence of the resummation in the infrared physics of 
bosons with the zero Matsubara frequency, while it vanishes just at the critical temperature Hr c = e T c 2 In 2 — /i B = 0. 
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E. Critical temperature 



The critical temperature in units of the Fermi energy tionship 2/x + = e T c 2 In 2: 
directly follows from Eq. (74) with the use of the rela- 



£ F 



2 
3^ 



2 7 -3-2m2 3C'(2) CW-2(ln2) 2 1 ^f^ 2 
8 2^2 3^2 + 8 U {~ 



= 0.260 - 0.0112 e + 0.0234 — + (9(e 2 ) 



In 2 e b 
3ir 2 £p 



(75) 



£ F 



where the numerical value Cn ~ 1.92181 is substituted. 
We find the critical temperature T c is an increasing func- 
tion of the binding energy near the unitarity limit. 
The next-to-leading-order correction is reasonably small 
compared to the leading term even at e = 1. The naive 
extrapolation of the critical temperature to the physi- 
cal case of d = 3 gives T C /£ F « 0.249 in the unitarity 
limit £b = 0. This value is surprisingly close to results 
from two Monte Carlo simulations, T c /ep — 0.23(2) [33] 
and T c /£p « 0.25 [36], while other two simulations pro- 
vide smaller values, T c /e F < 0.14 [34] and T c /sf = 
0.152(7) [35]. 

It is also interesting to compare the critical temper- 
ature in the unitarity limit with that in the BEC limit 
Tbec- I n the BEC limit, all fcrmion pairs are confined 
into tightly bound molecules and the system becomes a 
noninteracting Bose gas where the boson mass is 2m and 
the boson density is N/2. The critical temperature for 
the Bosc-Einstcin condensation of such an ideal Bose gas 
at d > 2 spatial dimensions becomes 



cbec 
£f 



r i + 



-2/d 



(76) 



To the leading and next-to-leading orders in e = 4— d, the 
ratio of the critical temperatures in the unitarity limit T c 
and in the BEC limit Tbec at the same fermion density 
is given by 



'- [1 + 0.0177 e + 0(e 2 



Tbec V 9 

= 0.943 + 0.0167 e+0(e 2 ). 



(77) 



The ratio is slightly below unity, indicating the lower crit- 
ical temperature in the unitarity limit T c < Tbec- The 
leading order term of the above ratio, T c /Tbec = \/8/9, 
has the following clear physical interpretation: The crit- 
ical temperature for the Bose-Einstein condensation at 
d = 4 is proportional to a square root of the boson's den- 
sity. In the BEC limit, all fermion pairs form the bound 
bosons, while only 8 of 9 fermion pairs form the bosons 
and 1 of 9 fermion pairs is dissociated in the unitarity 



limit [see the leading order terms of fermion and boson 
densities in Eqs. (66) and (69)]. Thus their ratio in the 
critical temperature should be T c /Tbec = \/8/9 < 1 at 
d = 4. 

The more appropriate estimate of T c at d — 3 will 
be obtained by matching the e expansion around four 
spatial dimensions with the exact result around d = 2. 
The critical temperature at unitarity in the expansion 
over e = d — 2 is given by T c = (e 7 /w) A, where A/e f = 
(2/e) e~ x l z is the energy gap of the fermion quasiparticlc 
at zero temperature [14]: 



T 

-L C 

£ F 



2eT" 



7T 



■ e- 1 '* [1 + 0(e)] . 



(78) 



We shall write the power series of e in the form of the 
Borcl transformation, 



T c (e) _ 2e<- x 



£F 



-V* / 
Jo 



dte-*B Tc {et). (79) 



Bt c (t) is the Borel transform of the power series in T c (e), 
whose Taylor coefficients at origin is given by Bx c (t) — 
1 + • • • . In order to perform the integration over t in 
Eq. (79), the analytic continuation of the Borcl transform 
B To (t) to the real positive axis of t is necessary. Here we 
employ the Pade approximant, where Bx c (t) is replaced 
by the following rational functions 



BtM 



A I 



1 + Pit j h PMt 

l + qit + --- + q N t N 



(80) 



Then we incorporate the results around four spatial di- 
mensions in Eq. (75) by imposing 



T c (2-e) 



£ F 



= 0.260-0.0112e+ ••■ 



(81) 



on the Pade approximants as a boundary condition. 
Since we have two known coefficients from the e expan- 
sion, the Pade approximants [M/N] satisfying M+N = 2 
are possible. Since we could not find a solution satisfying 
the boundary condition in Eq. (81) for [M/N] = [1/1], 
we adopt other two Pade approximants with [M/N] — 
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Tele, 




FIG. 7: (Color online) The critical temperature T c at uni- 
tarity as a function of the spatial dimension d. The upper 
solid line is from the expansion around d — 4 in Eq. (75), 
while the lower solid curve is from the expansion around 
d — 2 in Eq. (78). The middle two curves show the dif- 
ferent Borel-Pade approximations connecting the two expan- 
sions. The symbols at d — 3 indicate the results from 
the Monte Carlo simulations; T c /ef = 0.23(2) [33] (circle), 
T c /e F < 0.14 [34] (down arrow), T c /e F = 0.152(7) [35] (dia- 
mond), and T c /ef ~ 0.25 [36] (square). 



[2/0], [0/2], whose coefficients p m and q n are determined 
uniquely by the above conditions. 

Figure 7 shows the critical temperature T c in units of 
the Fermi energy e F as a function of the spatial dimen- 
sion d. The middle two curves show T c /ef in the dif- 
ferent Pade approximants connecting the two expansions 
around d = 4 and d = 2. These Borel-Pade approxi- 
mations give T c /e F = 0.173 and 0.192 at d = 3, which 
are located between the naive extrapolations to d = 3 
from the e = 4 — d expansion (T c /e F — > 0.249) and the 
e = d— 2 expansion (T c /ef — ► 0.153). It is also interesting 
to compare our results with those from the recent Monte 
Carlo simulations, where T c /e F — 0.23(2) [33], T c /e F < 
0.14 [34], T c /e f = 0.152(7) [35], and T c /e f « 0.25 [36]. 
Although these results from the Monte Carlo simulations 
seem not to be settled, the interpolation of the two expan- 
sions provides the moderate value T c /ef — 0.183 ±0.014 
not too far from the Monte Carlo simulations. 



F. Thermodynamic functions at T c 

Finally we show the thermodynamic functions at T c 
in the unitarity limit £b = to the leading and next-to- 
lcading orders in e. The pressure P normalized by the 
fermion density e F N follows from Eqs. (65), (73), and 
(74). Introducing the numerical values Cp « 8.4144 and 
Cn ~ 1.92181, we obtain the pressure up to the next-to- 



leading order in e as 
P 



£ F N 



0.116 + 0.0188e. 



(82) 



From the universal relationship in the unitarity limit E 
(d/2)P, the energy density is given by 



E 



e F N 



0.232 - 0.0205e. 



(83) 



The chemical potential at the critical temperature \x = 
eT c ln2 is O(e). Normalizing \x by the Fermi energy in 
Eq. (75), we have 



JL 



eln2 



2 
3^ 



0.180 e. 



(84) 



Then the entropy density T C S = (d/2 + 1)P— fJ>N is given 
by 



S_ 

N 



= 1.340- 0.642 e. 



(85) 



The next-to-leading-order corrections in the pressure and 
energy density are reasonably small compared to the 
leading order terms, while that is large for the entropy 
density. 

We match the thermodynamic functions at T c in the 
expansions over e = 4 — d with those around d = 2 as 
wc demonstrated for T c /e-p- The critical temperature 



around d = 2 is T c /ep 



-l/e 



which is exponentially 



small and negligible compared to any power series of 
e. Therefore, the pressure, the energy density, and the 
chemical potential at T c in the expansions over e = d — 2 
are simply given by those at zero temperature [14]: 



P 



e F N 



JL 

£T 



d + 2 



d + 2 



(86) 



(87) 



A straightforward calculation shows that the entropy per 
particle at T c to the leading order in e is given by 



S_ 

N 



3 £f 



2 7 re 7 - 1 



-X/e 



(89) 



Using the Borel-Pade approximations in order to con- 
nect the two expansions above, the thermodynamic func- 
tions at d = 3 are found to be P/(e F N)\ T — 0.172 ± 
0.022, E/(e F N)\ Tc = 0.270 ± 0.004, fi/e F \^ = 0.294 ± 
0.013, and S/N\^ = 0.642. The errors here indicates 
only the uncertainty due to the choice of different Pade 
approximants. In the recent Monte Carlo simulation, 
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the thermodynamic functions at the critical tempera- 
ture is P/(e F iV)| T = 0.207(7), E/(e F N)\ T = 0.31(1), 
H/e F \ Tc = 0.493(14), and S/N\ Tc = 0.16(2) [35]. We 
see that the interpolations of the two expansions indeed 
improve the series summations compared to the naive ex- 
trapolations from d = 4 or d = 2, while there still exist 
deviations between our results and the Monte Carlo sim- 
ulation. We can understand these deviations partially 
due to the difference in the determined critical temper- 
atures. The large deviations existing in /x/ep and S/N 
may be because we know only the leading term for /j, and 
the next-to-leading-order correction to S is sizable. 

V. SUMMARY AND CONCLUDING REMARKS 



In this paper, the thermodynamics of the Fermi gas 
near the unitarity limit at finite temperature has been in- 
vestigated using the systematic expansion over e = 4 — d. 
We discussed that the thermodynamics in the low tem- 
perature region T <C T c is dominated by the bosonic 
phonon excitations. The analytic formulas for the ther- 
modynamic functions at the fixed fermion density are de- 
rived in the two limiting cases; T <C eT c in Eqs. (41)-(44) 
and eT c «T« T c in Eqs. (51)-(54). 

In the high temperature region T ~ T c , the fermionic 
quasiparticles are excited as well as the bosonic quasipar- 
ticles. We showed that the similar power counting rule 
of e to that developed at zero temperature works even 
above T c . The superfmid phase transition at T = T c is 
of the second order. The critical temperature T c and the 
thermodynamic functions around T c were calculated to 
the leading and next-to-leading orders in e. We found 
the critical temperature is an increasing function of the 
binding energy £b near the unitarity limit: 

T £ 

— = 0.260 - 0.0112 e + 0.0234 — . (90) 

£F £F 

The next-to-leading-order correction is reasonably small 
compared to the leading term even at e = 1. In the 



unitarity limit £b — 0, the naive extrapolation of the 
critical temperature to the physical case of d = 3 gives 
T C /£ F « 0.249. 

We also discussed the matching of the e expansion 
around d — 4 with the expansion around d = 2. The 
critical temperature at unitarity in the expansion over 
e = d — 2 is given by T c /ef = (2e 7_1 /7r) e _1 / e , where 
its naive extrapolation to e = 1 gives T c /sf w 0.153. 
The Borel-Pade approximations connecting the two ex- 
pansions yielded T c /sf = 0.183 ±0.014 at d = 3, which is 
a moderate value located between the two naive extrap- 
olations (Fig. 7). These values are not too far from the 
results obtained by the recent Monte Carlo simulations 
where T c /e F = 0.15-0.25 [33, 35, 36]. We also applied 
the Borel-Pade approximations to the thermodynamic 
functions at T c , which yielded P/(e F N)\ T w 0.172, 
£/(£ F ^V)| Tc « 0.270, fi/e F \ Tc w 0.294, and" S/N\ n w 
0.642 at d = 3. 

The Borel-Pade approximations employed here to 
match the two expansions around four and two spatial 
dimensions do not correctly reflect large-order behaviors 
of the two expansions, i.e., the expansions over e = 4 — d 
and e = d— 2 are most probably not convergent and hence 
the Borel transform B{t) of such series expansions will 
have singularities somewhere in the complex i-plane [14] . 
In order for the accurate determination of the critical 
temperature and the thermodynamic functions at d = 3, 
it will be important to appropriately take into account 
the knowledge on the large-order behaviors of the expan- 
sions both around four and two spatial dimensions. The 
calculation of higher-order corrections to our results is 
also interesting for this purpose. These problems should 
be studied in future works. 
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